/*    Copyright 2010 Tobias Marschall
 *
 *    This file is part of MoSDi.
 *
 *    MoSDi is free software: you can redistribute it and/or modify
 *    it under the terms of the GNU General Public License as published by
 *    the Free Software Foundation, either version 3 of the License, or
 *    (at your option) any later version.
 *
 *    MoSDi is distributed in the hope that it will be useful,
 *    but WITHOUT ANY WARRANTY; without even the implied warranty of
 *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *    GNU General Public License for more details.
 *
 *    You should have received a copy of the GNU General Public License
 *    along with MoSDi.  If not, see <http://www.gnu.org/licenses/>.
 */

package mosdi.discovery;

import mosdi.distributions.PoissonDistribution;
import mosdi.util.Log;

/** Class that can precompute p-value bounds for Poisson distributions. */
public class PoissonBoundCalculator {
	private double minExpectation;
	private double maxExpectation;
	private double expectationRange;
	private int expectationSteps; 
	private double[][] table; 
	private boolean logarithmicTable;
	
	public PoissonBoundCalculator(double minExpectation, double maxExpectation, int maxValue) {
		this(minExpectation, maxExpectation, maxValue, 100, false);
	}
	
	public PoissonBoundCalculator(double minExpectation, double maxExpectation, int maxValue, int steps, boolean logarithmicTable) {
		Log.printf(Log.Level.DEBUG, "Pre-computing complementary Poisson CDFs in range %d to %d\n", 0, maxValue);
		Log.printf(Log.Level.DEBUG, "For %d different expectations from %f to %f, using %,d bytes\n",steps+1,minExpectation,maxExpectation, ((long)(steps+1))*(maxValue+1)*8);
		this.expectationSteps = steps;
		this.logarithmicTable = logarithmicTable;
		if (logarithmicTable) {
			this.minExpectation = Math.log(minExpectation);
			this.maxExpectation = Math.log(maxExpectation);
		} else {
			this.minExpectation = minExpectation;
			this.maxExpectation = maxExpectation;
		}
		this.expectationRange = this.maxExpectation - this.minExpectation;
		table = new double[expectationSteps+1][];
		for (int j=0; j<=expectationSteps; ++j) {
			double expectation = this.minExpectation+j*expectationRange/expectationSteps;
			if (logarithmicTable) expectation = Math.exp(expectation);
			table[j] = new PoissonDistribution(expectation).getCCDF(maxValue+1);
		}
	}

	/** Returns a value b, such that P(X>=k)>=b, where the distribution 
	 *  of X is Poisson with given expectation. */
	public double getLowerBound(double expectation, int k) {
		try {
			double x =  logarithmicTable?Math.log(expectation):expectation;
			int e = (int)Math.floor((x-minExpectation)/expectationRange*expectationSteps);
			return table[e][k];
		} catch(ArrayIndexOutOfBoundsException e) {
			throw new IllegalArgumentException("Argument not within range given on construction!");
		}
	}

	/** Returns a value b, such that P(X>=k)<=b, where the distribution 
	 *  of X is Poisson with given expectation. */
	public double getUpperBound(double expectation, int k) {
		try {
			double x =  logarithmicTable?Math.log(expectation):expectation;
			int e = (int)Math.ceil((x-minExpectation)/expectationRange*expectationSteps);
			return table[e][k];
		} catch(ArrayIndexOutOfBoundsException e) {
			throw new IllegalArgumentException("Given expectation not within bound given on construction!");
		}
	}

}
